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Abstract 

An essentially exact ground-state calculation algorithm for few-electron systems based on superposition of 
nonorthogonal Slater determinants (SDs) is described, and its convergence properties to ground states are 
examined. A linear combination of SDs is adopted as many-electron wave functions, and all one-electron wave 
functions are updated by employing linearly independent multiple correction vectors on the basis of the variational 
principle. The improvement of the convergence performance to the ground state given by the multi-direction 
search is shown through comparisons with the conventional steepest descent method. The accuracy and 
applicability of the proposed scheme are also demonstrated by calculations of the potential energy curves of few- 
electron molecular systems, compared with the conventional quantum chemistry calculation techniques. 
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Background 

In recent years, there have been many significant achieve- 
ments regarding electronic structure calculations in the 
fields of computational physics and chemistry. However, 
theoretical and methodological approaches for providing 
practical descriptions and tractable calculation schemes 
of the electron-electron correlation energy with con- 
tinuously controllable accuracy still remain as signifi- 
cant issues [1-15]. Although density functional theory 
(DFT) supplies a computationally economical and prac- 
tical method, there are many unexplored problems 
raised by unreliable results obtained for some systems 
in which highly accurate electron-electron correlation 
energy calculations are required, since results by DFT 
depend significantly on the exchange-correlation energy 
functional used to perform the calculation [16-18]. 

The available quantitatively reliable methods require 
higher computational costs than the DFT method [18]. 
Although quantum Monte Carlo methods [19-23] can 
be applied to molecular and crystal systems and show 
good quantitative reliability where extremely high- 
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accuracy calculations are required, difficulties in calcu- 
lating forces for optimizing atomic configurations are a 
considerable disadvantage and inhibit this method from 
becoming a standard molecular dynamics calculation 
technique. Configuration interaction (CI), coupled cluster, 
and Moller-Plesset second-order perturbation methods, 
each of which use a linear combination of orthogonalized 
Slater determinants (SDs) as many-electron wave func- 
tions, are standard computational techniques in quantum 
chemistry by which highly accurate results are obtained 
[24], despite suffering from basis set superposition and 
basis set incompleteness errors. The full CI calculation 
can perform an exact electron-electron correlation energy 
calculation in a space given by an arbitrary basis set. How- 
ever, it is only applicable for small molecules with modest 
basis sets since the required number of SDs grows explo- 
sively on the order of the factorial of the number of basis. 

The required number of SDs in order to determine 
ground-state energies can be drastically decreased by 
employing nonorthogonal SDs as a basis set. The reson- 
ating Hartree-Fock method proposed by Fukutome 
utilizes nonorthogonal SDs, and many noteworthy 
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results have been reported [25-30]. Also, Imada and co- 
workers [31-33] and Kojo and Hirose [34,35] employed 
nonorthogonal SDs in path integral renormalization group 
calculations. Goto and co-workers developed the direct 
energy minimization method using nonorthogonal SDs 
[36-39] based on the real-space finite-difference formalism 
[40,41]. In these previous studies, steepest descent direc- 
tions and acceleration parameters are calculated to update 
one-electron wave functions on the basis of a variational 
principle [25-30,36-39]. Although the steepest descent dir- 
ection guarantees a secure approach to the ground state, a 
more effective updating process might be performed in a 
multi-direction search. 

In the present study, a calculation algorithm showing 
an arbitrary set of linearly independent correction vec- 
tors is employed to optimize one-electron wave functions 
with Gaussian basis sets. Since the dimension of the 
search space depends on the number of linearly inde- 
pendent correction vectors, a sufficient number of cor- 
rection vectors ensure effective optimization, and the 
iterative updating of all the one-electron wave functions 
leads to smooth convergence to the ground states. The 
primary purpose of this article is to demonstrate the ad- 
vantage of using multiple correction vectors in searching 
for the ground state over the conventional steepest des- 
cent search in which only one correction vector is used. 
As a demonstration of the accuracy and applicability of 
the proposed calculation algorithm, essentially exact po- 
tential energy curves of few-electron molecular systems 
with long interatomic distances are described for cases 
where the conventional calculation methods of quantum 
chemistry fail. 

The organization of the article is as follows. In the 
'Optimization algorithm' section, the proposed calculation 
algorithm for constructing a basis set of nonorthogonal 
SDs by updating one-electron wave functions with 
multiple correction vectors is described. The expres- 
sion of the conventional steepest descent direction 
with a Gaussian basis set is also given for comparison. 
The convergence characteristics to the ground states 
of few-electron systems for calculations using single 
and multiple correction vectors are illustrated in 
the Applications to few-electron molecular systems' 
section. As demonstrations of the proposed calculation 
procedure, the convergence properties to the ground 
states of few-electron atomic and molecular systems 
are also shown. Finally, a summary of the present study 
is given in the 'Conclusions' section. 

Optimization algorithm 

The calculation procedures for constructing a basis set 
consisting of nonorthogonal SDs for AT-electron systems 
using single and multiple correction vectors are described 
here. An A/^-electron wave function \j/{ri, a^, i2' ff2.---> ^n> Ow) 



is expressed by a linear combination of nonorthogonal 
SDs as follows: 



\|/(ri,(Ti,r2,o-2, ...:rN,crN) 



^ CA<t^(ri, ffi, r2. 0-2...: Tn, (Tn) 



(1) 



Here, r, and Of denote the position and spin index of 
the i'th electron, respectively. L is the number of SDs, 
and ^'^(ri, ffj, r2, ff2,--; In< <^n) is the Ath SD, given by 
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(3) 



with (pf{i) and fi(cr,) being nonorthogonal and 
unnormalized one-electron basis functions and spin 
orbital functions, respectively. The one-electron wave 
function (pf(i) is constructed as a linear combination 
of Gaussian basis functions Xs(t) [24] as 



(4) 



Here, M and Df^ are the number of basis functions 
and the sth expansion coefficient for the rth one- 
electron wave function 0f (r), respectively. 

The steepest direction is implemented in the expression 
of the total energy functional E of the target system on the 
basis of the variational principle, without the constraints 
of orthogonality and normalization on the one-electron 
wave functions. The updating procedure of the pth one- 
electron wave function belongs to the Ath SD which is 
represented as 



ryt (new) _ (old) , A j^A 



(5) 



where cC^ is the acceleration parameter, which is deter- 
mined by the variational principle with respect to the total 
energy E, i.e., [28] 



dE 
da^ 



(6) 
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The component of the steepest descent vector Kp^y„ is 
given by 



where 



^J^-lS^B^fs^'-'TD^ (y \y)S 



" p,m J I A=\B=\ 



(7) 



(8) 



where the elements Sf^ are the overlap integrals between 
the one-electron basis functions, i.e., 



Sf = \dr cpf *(r) <(r). 



(12) 



{H^) and {Hf^) are the matrix elements of the 
Hamiltonians, 



^o = E 
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(13) 
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Here, denotes the element of the /th row and /th 

column of the matrix S^^ . The overlap integral matrix 
S is defined by 



rAB cAB oAB 
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and 



N N 

■^^ = 2^^lr -r I' 



(14) 



respectively. Here, V(r) stands for an external potential. 

The proposed calculation procedure employs 
linearly independent multiple correction vectors for 
updating the one-electron wave function. The jCth 
one-electron wave function in the ^th SD is updated 
by 



pA (new) _ (old) , ^ r,A 



(15) 



where Cj(j = 1, 2,..., L + N^) and are the expansion 
coefficient and the number of correction vectors, re- 
spectively. The components of the correction vectors 
G^,m determine linearly independent correction 
functions ^^(r) which are defined as linear combinations 
of Gaussian basis functions as 



(16) 



Since the linearly independent correction vectors can be 
given arbitrarily, randomly chosen values are employed in 
the present study. A larger number of correction vectors 
Nc realize a larger volume search space; however, the 
number of the linearly independent vectors Nc is restricted 
to the dimension of the space defined by the basis set 
used. 
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Thus, we have a linear combination of L + Nc SDs as 
the new A/'-electron wave function 



\|/(ri,o-i,r2,o-2, ■■■,rM,(rN 
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where 



L+Nc 

(17) 



(18) 



Figure 1 illustrates the flow of the present calculation 
procedure. Unrestricted Hartree-Fock (UHF) solutions 
for a target system are used for initial one-electron wave 
functions. The coefficients CjlJ = 1,2, ...,L + N^) of 
Equation 17 are given by solving the generalized eigen- 
value equations obtained by employing the variational 
principle applied to the total energy, and we can have a 
new TV-electron wave function as a linear combination of 
L SDs as shown in Equation 17. Iteration of the above 
updating process for all the one-electron wave functions 
of all SDs increasing the number of the SDs' L leads to 



[ Hartree-Fock (L=l) ] 



P=1,...,N 
A=1,...,L 



Matrix elements I'i^ 



H 



Eigenvalue eq. HC=ffSC 




[ END I 
Figure 1 Flow of the present algorithm. 



an essentially exact numerical solution of the ground 
state. 

Applications to few-electron molecular systems 

Convergence performances for searching for the ground 
state of a C atom with the 6-31G** basis set are shown 
in Figure 2. The UHF solutions are adopted as initial 
states, and the number of employed SDs is 30. The 
steepest descent direction and acceleration parameter 
are adopted for the calculation using one correction vec- 
tor (Nc =1), and seven randomly chosen linearly inde- 
pendent correction vectors are added to the steepest 
descent correction to create a calculation with eight 
correction vectors (Nc =8). An indispensable advantage 
of the multi-direction search over the single steepest 
descent direction search is clearly demonstrated. Al- 
though the steepest descent vector gives the direction 
with the largest gradient, it does not necessarily point 
toward the global energy minimum state. On the con- 
trary, a linear combination of multiple correction vectors 
can be used to obtain the minimum energy state within 
the given space by adopting the variation principle. 

Figure 3 illustrates the convergence performance of the 
proposed method for the electron-electron correlation 
energy of a HF molecule with the 6-31G** basis set as a 
function of the number of employed SDs. Calculated 
correlation energies are shown by ratios to exact ones 
obtained by full CI. The convergence performance to 
the exact ground state is improved by increasing the 
number of correction vectors, since the volume of the 
search space for a one-electron wave function increase 
with increasing Nc- The essentially exact ground-state 
energy is obtained using less than 100 nonorthogonal 
SDs with an error of 0.001%, compared with the exact 
value in which 99.5% of the electron-electron correl- 
ation energy is counted. The obtained convergence is so 





-37.67 n 




-37.68 , 




-37.69 - 


3 


-37.70 - 


(d 


-37.71 - 




-37.72 - 


<D 
C 


-37.73 - 


LU 


-37.74 - 




-37.75 - 




-37.76 - 



C 6-3 1G 



/Vc=1 
Nc=B 



Exact 



)-o-o-o-o-o-o-o-o 



0 5 10 15 20 25 30 
Number of iterations 

Figure 2 Effectiveness of multi-direction search on total energy 
convergence. Effectiveness of multi-direction search on total energy 
convergence as a function of the number of iterations for a C atom 
with the 6-31G** basis set is shown. 
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Figure 3 Convergence performance of the proposed method 
for the correlation energy. Convergence performance of the 
proposed method for the correlation energy of a HF molecule with 
the 6-3 IG** basis set as a function of the number of employed SDs 
is shown. 



smooth that the accuracy of the total energy is control- 
lable by adjusting the number of employed SDs. On the 
other hand, the full CI method requires over 10^ orthog- 
onal SDs, and thus the reduction in the numbers of SDs is 
a significant advantage of adopting nonorthogonal SDs. 
The ground-state energy obtained by the proposed 
method does not depend on the components of the 
correction vectors; however, the rate of convergence 
does depend on the number of employed correction 
vectors N^. 

The potential energy curve calculated when a single 
H atom is extracted from a CH4 molecule as shown in 
Figure 4. Calculations are performed using the 6-31G* 
basis set. Although the bond lengths are close to the 
equilibrium one, the errors in the energies obtained by 
coupled-cluster theory with singles and doubles (CCSD) 
plus perturbative triples (CCSD(T)) are a few milliHartree; 
at longer bond lengths, the accuracy of the results 
appears to deteriorate [42]. In contrast, the proposed 



calculation procedure ensures essentially exact ground 
states at all bond lengths, since no approximations are 
employed. 

Figure 5 illustrates the potential energy curve along 
the symmetric stretching coordinate of a H2O molecule 
in the 3-21G basis set. The angle between the O-H 
bonds is fixed at 107.6°. These results shown for the pro- 
posed calculation method, CCSD and CCSD(T) exhibit 
the same trends as for a CH4 molecule. The results for 
near the equilibrium bond length demonstrate comparable 
performance between the four methods, whereas results 
for long bond lengths indicate only that the proposed 
method has comparable performance with full CI not 
producing the same unphysical energy curves as CCSD 
and CCSD(T) around 2.3 A [42]. 

Conclusions 

A reliable and tractable technique for constructing the 
ground-state wave function by the superposition of 
nonorthogonal SDs is described. Linear independent 
multiple correction vectors are employed in order to up- 
date one-electron wave functions, and a conventional 
steepest descent method is also performed as a comparison. 
The dependence of convergence performance on the num- 
ber of adopted correction vectors is also illustrated. The 
electron-electron correlation energy converges rapidly and 
smoothly to the ground state through the multi-direction 
search, and an essentially exact ground-state energy is 
obtained with drastically fewer SDs (less than 100 SDs 
in the present study) compared with the number re- 
quired in the full CI method. For the few-electron 
molecular systems considered in the present study, 
essentially exact electron-electron correlation energies 
can be calculated even at long bond lengths for which 
the standard single-reference CCSD and CCSD(T) show 
poor results, and the practicality and applicability of the 
proposed calculation procedure have been clearly dem- 
onstrated. In future studies, calculations employing 
periodic boundary conditions and effective core potentials 
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Figure 4 Potential energy curve of a CH4 molecule obtained 
using the proposed algorithm with 6-3 1G* basis set. 
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Figure 5 Potential energy curve of a H2O molecule obtained 
using the proposed algorithm with 3-21 G basis set. 
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(ECPs) [43] will be performed. A new procedure to reduce 
the iteration cost should be found in order to increase the 
applicability of the proposed algorithm for the calculation 
of essentially exact ground-state energies of many-electron 
systems. 
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